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We investigate the effect of phenomenological relaxation parameters on the third order optical 
nonlinearity of doped graphene by perturbatively solving the semiconductor Bloch equation around 
the Dirac points. An analytic expression for the nonlinear conductivity at zero temperature is 
obtained under the linear dispersion approximation. With this analytic formula as starting point, 
we construct the conductivity at finite temperature and study the optical response to a laser pulse 
of finite duration. We illustrate the dependence of several nonlinear optical effects, such as third 
harmonic generation, Kerr effects and two photon absorption, parametric frequency conversion, 
and two color coherent current injection, on the relaxation parameters, temperature, and pulse 
duration. In the special case where one of the electric fields is taken as a dc field, we investigate 
the dc-current and dc-field induced second order nonlinearities, including dc-current induced second 
harmonic generation and difference frequency generation. 
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I. INTRODUCTION 

The optical nonlinearities of monolayer graphene have 
recently attracted wide attention/^— both experimen¬ 
tally and theoretically. The nonlinear susceptibility of 
graphene^ is both strong - per atom it is orders of 
magnitude higher than that of common gapped semi¬ 
conductors and metals - and controllable by the chem¬ 
ical potentialf^i^ which can be tuned by an external 
gate voltage^i^ or chemical dopingi^ With the possibil¬ 
ities it offers for integration in silicon-based optical in¬ 
tegrated circuits, graphene is an exciting new candidate 
for enhancing nonlinear optical functionalities in silicon- 
based on-chip optical devices, such as on-chip broad¬ 
band light sources, electro-optic modulators /°d^ opti¬ 
cal switches/^4ii and optical transistorsIn realiz¬ 
ing some of these devices/^ the presence of second order 
optical nonlinearities, especially second harmonic gener¬ 
ation (SHG), is a key requirement. 

The third order optical nonlinearity is described by the 
susceptibility tensor W2, W3) or equivalently the 

conductivity tensor (wi, a;2, W3), which has a complex 
frequency dependence. It describes different physical ef¬ 
fects, such as third harmonic generation (THG), which 
is determined by w, w); Kerr effects and two pho¬ 
ton absorption, which are determined by w, w, w); 

two-color coherent current injection, which is deter¬ 
mined by w, —w, 2a;); and parametric frequency 

conversion (four wave mixing), which is determined by 
—u}s,u}p,ijjp). Due to the inversion symmetry of 
its crystal structure, pristine graphene has no second 
order optical nonlinearities arising from electric dipole 
transitions. However, in graphene-based photonic de¬ 
vices an effective second order susceptibility can arise 
from the breaking of inversion symmetry in a number 


of ways: (1) the presence of an asymmetric interface be¬ 
tween graphene and the substrate/^^S^ not relevant for 
normally incident light; (2) the contribution of forbidden 
transitions involving the finite wave vector of light 
(3) the presence of natural curvature fluctuations of sus¬ 
pended graphene^ (4) the application of a dc electric 
field to generate an asymmetric steady statei^d^ ~^^:^^i^^ 
The last is associated with the third order optical non¬ 
linearity x*'^H'^ 1 i‘^ 2,0), with one of the electric fields in¬ 
dependent of time. It includes current induced second 
harmonic generations^ (GSHG) or electric field induced 
second harmonic generation (EFISH). 

Experimental studies of many of the optical nonlinear 
effects mentioned above have already demonstrated in 
graphene. Typically the experimental data are analyzed 
by extracting an effective optical nonlinear susceptibil¬ 
ity, with the graphene monolayer treated as a thin film 
with a thickness of 3.3 In this way, most of the 

experimental techniques used to determine the nonlinear 
optical response of bulk materials or thin films can be 
directly applied to the study of graphene. In a gapped 
semiconductor, third order susceptibilities do not change 
drastically in the nonresonant regime, where all photon 
energies are much lower than the energy gapj^ Yet they 
show a strong and complicated photon energy depen¬ 
dence in pristine graphene because resonant transitions 
always exist for any photon energy, due to the vanish¬ 
ing gap and the presence of free carriers, leading to some 
similarities with a metal filmi^ These complexities have 
been observed in experimental studies of parametric fre¬ 
quency conversioufSi THG /^1^^1^^ Kerr effects and two 
photon absorption)^i2&4^ two color coherent control^^^li 
and in graphene. 

Theoretically, the optical nonlinearities of graphene 
have been investigated by perturbative treatments based 
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on Fermi’s Golden Rule, and by density matrix calcu¬ 
lations, both of which are standard methods in study¬ 
ing the optical response of gapped semiconductors. In 
an earlier communication we sketched some of the rele¬ 
vant work done before early 2014^; recent contributions 
include a calculation by Mikhailov^ of THG^, and nu¬ 
merically solution of the equation of motion under strong 
laser fields by Avetissian et All of these studies 

focused on one or a few nonlinear effects. In our earlier 
worb^ we performed a perturbative calculation based on 
a density matrix formalism; ignoring all scattering ef¬ 
fects, we obtained an analytic expression for the general 
optical sheet conductivity which can be 

related to the effective susceptibility ^2,^3), in 

doped graphene at zero temperature. We found that the 
optical conductivities depend strongly on the chemical 
potential and photon energies, and exhibit many diver¬ 
gences associated with resonant transitions, which occur 
when photon energies or their combinations match the 
chemical potential gap. Taking W3 = 0 and including 
phenomenological relaxation times for the generation of 
both dc and optically induced current, we calculated the 
current induced second order nonlinearities at zero tem¬ 
perature and obtained an analytic expressioi>S for CSHG. 
The effective susceptibility shows two peaks correspond¬ 
ing to two resonant transitions induced by the fundamen¬ 
tal and the second harmonic light, with the peak values 
strongly dependent on the relaxation time. Adopting the 
parameters used in calculations of bilayer graphene?^ we 
obtained a prediction of a peak susceptibility in mono- 
layer graphene that was similar to that predicted for the 
bilayer; the EFISH contribution was ignored in that cal¬ 
culation. 

The importance of the relaxation time demonstrated 
in that study, and the desire for more realistic calcula¬ 
tions to compare with experiment, motivates the present 
work. Here we consider the inclusion of scattering effects 
in the semiconductor Bloch equations (SBE) within a re¬ 
laxation time approximation, allow for finite temperature 
to the extent that it affects the initial state, and explic¬ 
itly consider the nonlinear response to pulses of light. We 
obtain an analytic expression for the full nonlinear op¬ 
tical conductivity W3) for optical transitions 

around the Dirac points. We discuss the predictions that 
follow from this expression for different optical effects, 
and we compare with experiment where possible. 

Our focus in this work is on doped graphene, where 
the chemical potential ^ 7^ 0. However, the chemi¬ 
cal potential dependence of our general expression for 
allows us to study the special case of 
/i —0. At the very least we might expect that, for 
electrons close to the Dirac points, the distinction be¬ 
tween “interband” and “intraband” motion could be lost. 
Although different terms that are nominally associated 
with interband and intraband motion arise naturally in 
the development of the perturbation series, the distinc¬ 
tion between those two “kinds” of motion is at best 
approximate;^ and we indeed hnd that the way those dif¬ 


ferent formal terms contribute to the final result for small 
fj, is nontrivial. More importantly, we generally associate 
the validity of a perturbative expansion of the optical re¬ 
sponse with the assumption that the energy induced by 
the presence of the optical field is much less than the en¬ 
ergy difference between the bands. In graphene this is 
always violated for some states around the Dirac points, 
regardless of the strength of the optical held. If these 
states are occupied by electrons, as they are in undoped 
graphene, the reasonableness of a perturbative expansion 
is in doubt. Indeed, even a semiclassical treatment of the 
response to an applied electric held of electrons near the 
Dirac points exhibits a breakdown of the perturbative 
analysis^^ as /r —>■ 0. We hnd evidence of the same kind of 
behavior in the quantum treatment presented here. This 
has consequences even for doped graphene if hnite tem¬ 
perature is considered, for thermal huctuations always 
place some electrons near the Dirac points. 

We organize our paper as follows: In Section |H] we 
introduce the SBE and our approximations for includ¬ 
ing scattering effects; the details of the derivation of the 
nonlinear optical conductivity is given in Appendix 
The last two subsections of Section El address the ex¬ 
tension of the calculation to hnite temperature, and the 
treatment of the response to a pulse with hnite dura¬ 
tion. In Section uni we discuss the third order nonlin¬ 
ear effects, including THG, Kerr effects and two pho¬ 
ton absorption, two-color coherent current injection, and 
parametric frequency conversion; in Section IIVI we dis¬ 
cuss the current-induced second order nonlinearities, in¬ 
cluding CSHG, EFISH, and the nonlinear optical con¬ 
ductivity a;s,a;p,0) that describes current-induced 

diherence frequency generation. Throughout the sections 
we compare with experimental results when appropriate. 
We conclude in Section 0 

II. MODEL 

We take the Hamiltonian of graphene to be 

H = Hq + H(.fi + Hgp + Hgi + I (1) 

Here Hq is the unperturbed electron Hamiltonian, 

^0 — ^ ^ / dkc^y,(xsii , (2) 
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where the ask are annihilation operators of Bloch states 
\sk) for band s and wave vector fc, with eigen energy 
Esk- Here H^r describes the interaction with radiation 
and in the dipole limit, where the electric held E{t) is 
approximated as uniform, we have 

dleR — eE/(t) ■ ^ ^ / dka^^ ^ (^si ^ k^ ^Sik : 

S 1 S 2 

(3) 

where e = — |e| and 

f dv 

6iS2fe = M — Us^kir)V kUs^k{r) (4) 

^cell ‘^cell 
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is the Berry connection between states |sifc) and |s2fc), 
with ^ceii the unit cell area and Usfc(r) the periodic 
part of the Bloch function, (r|sfc) = ( 27 r)“^e®^''^Msfc(r, z), 
where k = {kx,ky) and r = {x,y); the graphene is as¬ 
sumed to lie in the x — y plane. We neglect any response 
of the system to electric field components in the z di¬ 
rection. The scattering terms are given by Hgi for the 
electron-impurity scattering, H^p for the electron-phonon 
interaction, and Hge for the carrier-carrier scattering. 

The system is described by a density matrix that is 
initially diagonal both in band index and (continuous) 

wave vector, = nsikiSs:^s 25 {ki - k2), 

where 0 < ^ 1 describes the initial occupation of 

the state. In the presence of an applied uniform electric 
field it remains diagonal in fci and k2 but can acquire 
off-diagonal elements in si and S2-, describing the cor¬ 
relation between state amplitudes for |sifei) and 152^2), 

= Psis2kiit)S{ki - k2). We can think of 

PsiS2k{t) as the elements of a 2 x 2 matrix pk(t), and their 
dynamics are determined by the SBE 


^ dPsiS2k ^ _ es2k)PsiS2k 


dt 


-\- ieE (t) ^ ( {^siskPsS2k PsiskCsS2k) 


- eE^{t) 


9 psiS 2 k 

Oka 


9 psiS 2 k 


dt 


( 5 ) 


Here 


dPsi S2fc 
dt 


scat 


includes the scattering terms induced by 


Hei+Elep+Heei which could in principle be obtained from 
well-established treatments of many-particle systems, 
such as the many-particle density-matrix frameworb^Si^ 
or the Keldysh Green function method;?^ In an ordi¬ 
nary semiconductor with parabolic band structure, the 
current relaxation is mostly caused by carrier-phonon 
and carrier-impurity scattering, while carrier-carrier in¬ 
teractions are less significant due to the approximate 
equivalence of momentum conservation and velocity con¬ 
servation. However, the novel linear band structure 
of graphene breaks this equivalence, and the carrier- 
carrier interactions play an important role in current 
relaxationj^ii^ thus the full expression for the scat¬ 
tering terms is complicated and even hard to solve 
numerically^ 

We proceed in the standard way by assuming the va¬ 
lidity of a perturbation expansion 


ps^s 2 k{t )= 


( 6 ) 


n=0 


oc E". Here = 5 s,s 2 ns,k 

is the density operator characterizing the equilibrium oc¬ 
cupation of single-particle states at finite temperature T 
and chemical potential p, Usk = [1 + is 

the Fermi-Dirac distribution with /3 = 1 /{ksT) where ks 


is Boltzmann’s constant. From Eq. satisfies 


t- ^PsiS2fc _ ./ \ (n) 

n — t{Ssik £s2fc)PsiS2fe 


dt 


zeE<^{t)J 2 {i 

S 

S 2 k 


(n—1) 
siskP ss 2 k 


rs\sk SsS2/c J 


- eE^{t) 




dp] 


(n) 

S1S2/C 


dka 


dt 


( 7 ) 


where = 0 for n < 0 . As a very rough approximation, 
a relaxation time approximation^ can be adopted to give 


dpl,S 2 k 


dt 


_ _p(") (n) 

^ SiS2krsiS2k ’ 


for n > 1. 


( 8 ) 


(n) 

Here is a relaxation parameter introduced to de¬ 
scribe the dynamics of corre¬ 

sponds to a phenomenological relaxation time. In a real 
system, can be expected to depend on the tem¬ 

perature, chemical potential, and external field.-^ Yet 
because the relaxation plays an important role in optical 
nonlinearities around resonant transitions, the extremely 
phenomenological treatment^ in Eq. can still reveal 
part of the physics, and in a very simple way. Even with 
the use of the six phenomenological constants T^”^ = ri”^ 
for interband transitions and T^”], = T^”^ for intraband 
transitions, we are still able to obtain an analytic result 
for the perturbation calculation within the linear disper¬ 
sion approximations around the Dirac points at zero tem¬ 
perature. From p^i^\t), the (areal) current density, which 
in our model has only x and y components, is calculated 
as J^(t) = with 


j(");'^(t) = e ^ 


( 9 ) 


We give the derivation in Appendix where the spin 
degeneracy is included. We extract the linear optical 
conductivity cr^^^’‘^“(a;) from 

j(i);d(i) = f J^CT(i):‘'“(w)E“(a;)e-*‘"*, (10) 

J ZTT 


where E°'{ui) = f dtE“(t)e“‘. In graphene, the hexag¬ 
onal lattice has Eek (6/mmm) symmetry;^ and there 
is only one independent nonzero component = 

fjd)',yy^ We first consider the zero temperature results. 
In this paper, we restrict ourselves to the neighborhood 
of the Dirac points (see Fig.®, assuming a linear disper¬ 
sion relation with two relevant bands that we label s = -I- 
(upper) and — (lower). We recover the usual result^i^ 




^\p\ 

fvjj + frf) 


( 11 ) 
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FIG. 1. (color online) Illustration of the linear dispersion ap¬ 
proximation of the graphene band structure around the Dirac 
point. The arrows show optical transitions induced by one 
photon with energy hw (right) or three photons with energy 
hwi (left). 

Here cto = e^/( 4 /i) is the universal conductivity, and 
with 1? = iJr -b i'di is 


(s') 

d = fvjj + iTe \ Wo = W2 +W3, and w = wi -bwo. We have 
followed the standard convention of nonlinear optics^ in 
symmetrizing the terms (w^,Wfe,w/) by permut¬ 

ing the indices {jkl) to arrive at the nonlinear conduc¬ 
tivity cr(^)’‘^“^‘^(wi,W2,W3). The light-matter interaction 
in Eq. ([ 3 ]) can be formally separated into an interband 
contribution (si yb §2) and an intraband contribution 
(si = 32), and the terms proportional to the different St 
in CT(3);dabc classified according to how many times 

each contribution appears 4 S The term proportional to 
iSi arises from only the intraband contributions, and the 
term proportional to Sq arises from only the interband 
contributions; all others involve mixtures of both. The 
quantities -^{ 3 )-,dabc^ ^{ 3 )-dabc gj.g fourth or¬ 

der tensors. Neglecting the optical response in the z di¬ 
rection, there are in all 8 nonzero components for the 
symmetry, among which three are independent; they are 

fj{3y,xxyy _ ^{3)\yyxx 
fji3)-,xyxy _ ^{3)\yxyx 
^{3)-,xyyx _ ^{3)-,yxxy ^ 


= In 


2|/r| -|- "d 


2 \fj.\ -19 

'dr “b 2|^| 


+ i { TT + arctan 


dr - 2|/r| 
■di 


— arctan 


d,. 


( 12 ) 


The third order current is given as 

X E;“(wi)E''(w 2 )E'=(w 3 )e-*(‘^i+“=+‘^ 3 )‘, ( 13 ) 


Here the symmetrized third order optical conductivity 

^{3);dabc -g 


^( 3 );dabc(^^,^ 2 ,W 3 ) 

1 -b ct(^)’‘'^“(w 2, W3, Wi) 

6 

+5:(^);‘^“*'(w 3, Wi, W2) -b W3, W2) 

^~( 3 );dc 6 a(g^ 3 ^ t^2, Wi) + ^ 


,(14) 


where the unsymmetrized third order optical conductiv¬ 
ity is given as 


^dabc 


= 


+ 


^dabc 


ids) , St^^ido) 


+ 


+ 


VVqV^ Wq 

St^^{do,d 3 ) Si^^%d) St^^{d,d 3 ) 


^dabc 


1 / 

{d,do) 


:'dabc / 


^0 


V3 


+ St'>^{d,do,d 3 ) 


( 15 ) 


with (T3 = cro(/iUFe)^/7r, z/3 = SW3 -b fT^, -^3 = Sws -b 
fq = hjjJ3 + iTf \ do = huJo + ir''e\ V = hjjj + iTf\ 


and 


^(3)-,xxxx _ ^(3)-,yyyy 


j3)-,xxyy , (3)-,xyxy 


+ a 


{3)-,xyyx 


( 17 ) 


In the following, we write the independent nonzero com¬ 
ponents of fourth rank tensors as column vectors, order¬ 
ing the independent components of a fourth rank tensor 

'rp{3)-,xxyy~ 

By employing the constant 


Tdabc gg y = 

vectors 


rp{3)-,xyxy 

J^i3)-,xyyx 



-3 


1 


■ 1 ■ 


'l' 

II 

1 

1 

to 

II 

-3 

1 

II 

CO 

1 

-3 

II 

0 

1 

1 


( 18 ) 

where note Aq = —(Hi 4-^24-^3), we can present the an¬ 
alytic expression for the different components of '^A)Aa-bc 
appearing in Eq. (flSl) at zero temperature, using the ap¬ 
proximation of a linear dispersion relation around the 
Dirac points, as 


' 

H3 1 H3 


S 3 {do) = 'H^{do)^-—^ 
Siido.do) = -gy^ido) 


do ImI do 

d 3 A 2 4 - 'd2H3 


dldl 


-Uyido) 


{do + d 2 )A 2 + d 2 A 3 


A2 


dod 2 


dld ^2 

1 ^3 

ImI dodo 


( 19 ) 

( 20 ) 
( 21 ) 


( 22 ) 
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( 23 ) 


S,{^,^ 3 ) = -G,{^z):^T 7 ^ 


+Gf^W 


i^o) = 'Hfj.i'i^o) 


^2(^2 -^2) 

^ 3^0 , 1 Aq 


— d1) \fi\i)'d 3 


( 24 ) 


Ao 


^3 


+ K.W I # - ^ 


Ao 1 -^3 


(25) 


and 


Ssi'd,'do,i^ 3 ) 




[(^-^ 3 )^ 2^3 

+ 1^3 + l?l? 3 (- 3 l ?0 + 2 i? 3 ) 


+ Gfj.i'&o) 


G>.m 


'do'&l + '&1'&2 — '&0'&2 . 

^PPI ^ 

dl'&l'&l 

1 


i 9 i?f(d-)d 3 ) 

+ ^3(^0 + T^s) - ^?( 3 i?o + 4 i? 3 ) 


+ -Hf^i^) 


Ao 


_ ^3)3 

2^3 


^43 


4d^ — 3i?t?o — 2^?■d3 + i?o'!^3 


i 9 i?f (i? - Z? 3 )^ 


213 


+ M^) 


1 243 


|/z| ' 

where i?2 = 1^0 — i^3j i?i = 'd — 'doi and 


= 


2|/r| — i? 2\fi\ + i 9 ’ 

1 1 


(26) 

( 27 ) 

( 28 ) 


( 2 | m |+^)2 ( 2 | m |-^)2 

For the details see Appendix El 

Using the nonzero independent components, the third 
order current in Eq. (d can be written as 

(3) _ f dwidw2dw3 

^’~J ( 27 r )3 

X [a(3);""^^(wi,W2,W3)-E(wi)£;(w2) • Eitos) 

+a^^y’^y^v(u;3,u;2,u;3)E{0J2)Eiu;i) ■ Eicos) 
+a^^)-^yy^{u;i,u} 2 ,u; 3 )E{oj 3 )EiiJi) ■ E{u;2 )] • ( 29 ) 


A. Divergences and limits 


The results for d^^y’'^°-’^‘^(uji,uj2,0J3) show a complicated 
dependence on the Wj, on the and on /i. The ex¬ 
pressions in Eqs. (jl 9 ll 26 l) seem to exhibit a number of 
divergences, but some of them are only apparent: For 
example, there seem to be divergences when d — -03 = 0, 
but a careful collection of terms shows that even in 
the absence of relaxation lim5_>.o (t( 3 );‘^“^'=(—w, w -|- 6,^3) 
is finite. Some of the divergences are of course real: 
There are divergences for 2 |/r| ± = 0 in the func¬ 

tions G (d), % (d), and X (d), which lead to divergences 
in cf^3)Ra6c(^3^^2,a;3). These are associated with inter¬ 
band optical transitions, and for nonvanishing relaxation 
they occur at frequencies removed from the real axis; 
we will see how some of them affect the structure of 
cr(3);‘^“^^(a;i,a;2,a;3) in Sections Hill and ITVl There are also 
divergences associated with d -|- ds = 0 . In the absence 
of relaxation these occur when wi -I- a;2 -I- 2023 = 0 , and 
lead to a divergence in (T(3);da6c^_^^ _ 2^3 _|_ J^a;2,a;3) 
as A special case of these is when d = 0 and 

dj = 0 for j = 0 , 1 , 2 , or 3 . Some of the associ¬ 
ated conductivity terms, such as w, w, w) and 

cj(dydabc^_^^ —a;, 2 a;) will be considered in Section uni 
All of these divergences only occur at complex frequen¬ 
cies in the presence of relaxation, and have their analogs 
in gapped systems. Of a different nature are the diver¬ 
gences that arise as |/i| —>■ 0 . While in a semiclassical cal¬ 
culation and in the absence of relaxation the intraband 
third order nonlinear response coefficient that can be ex¬ 
tracted from the full nonlinear response is divergent!^ as 
\^\ , one might hope that in the presence of relaxation 

this would be ameliorated. Yet in general it is not. To 
see this, we reorganize the unsymmetrized conductivity 
to write 

CT(^);‘^“''^(a;i,a;2,a;3) = a;2, a;3) 

, •~( 3 );(ia 6 c / \ 

+ 0-3 {Wi,W 2 ,W 3 ) 

+ 'd^c’’d°‘^''{uji,uj2,i03) ( 30 ) 

where a a includes all terms involving t/^, ds includes 
all terms involving and and the remainder, ctc, 
includes all terms proportional to |/r|“ 3 . Similar sepa¬ 
rations are also used for the symmetrized conductivity 
^(3);dabc^ The term cfq can be simplified to yield 


~(3)-,dabc, X 0-3 /Aq A 3 \ 

vdv3d3 


rf^) 


.( 31 ) 


Note that even for finite relaxation we have 
o'^)’‘^“^^(aJi,a;2,a;3) diverging as |^| —>■ 0, for gen¬ 


eral frequencies (a;i,a;2,a;3), when ^ F)'^"' for both 
j = 1 and j = 3 . At least within the simple description of 
relaxation we adopt here, the perturbation theory seems 
to demand that either the first or third order relaxation 


dj) 
































6 


rates (or both) must not distinguish between intraband 
and interband relaxation to achieve a finite result as 
1^1 —>■ 0 . This is at least consistent with the physical 
intuition that the distinction between intraband and 
interband motion is blurred as |/r| —>■ 0, in any case for 
electrons near the Fermi level, and any reasonable theory 
should respect that; recall that in our phenomenological 

description of relaxation all carriers share the same F^'^^ 
(i) 

and F) . But clearly a more sophisticated theory is in 
order to address the limit |^| —>■ 0. 

More evidence for the blurring of the distinction be¬ 
tween intraband and interband motion as \fi\ —>■ 0 can be 
seen from how the contributions to W3) 

arise. The term in a;2,wa) that contains 

only contributions from the formal intraband (si = S2) 
component of Eq. (jS]) is the term proportional to 5 i; it 
varies with |/i| as |^| , which is qualitatively different 

than the variation as |/r| of the corresponding Drude term 
in the linear conductivity. Yet as |p,| 0 the contri¬ 
bution to clI 2,wa) involving only the formal 

interband (si yf S2) component of Eq. ([ 3 |), that is pro¬ 
portional to iSe, also becomes important; while it includes 
contributions from C/^(0a) and Q^{0), there is also a term 
proportional to |/i| The formally “mixed” terms, Sj^ 
with j = 2 , 3 , 4 , 5 , 7 , 8, also provide terms proportional 
to The summation of all these terms, all formally 

involving different proportions of interband and intra¬ 
band contributions, gives Eq. dSID; the 1^1 ^ behavior in 
(Tc (wi, a;2, W2) cannot be associated with motion that is 
just formally intraband. 

Now note that crc(wi, W2, W2) vanishes as all the F^"g 
vanish. Yet here we physically would expect to recover 
the relaxation free, semiclassical result^^ of a perturba¬ 
tive response divergent as |/r| for huji <C |p.|, and the 
result associated in that calculation with purely intra¬ 
band motion. And we do recover it here, but in a non- 
trivial way: Although {uji,uj2,uj3) vanishes, when the 
other contributions to a;2, wa) are assembled 

and the limit hiOi <C |^| taken we find 


a;2, wa) the results derived earlier— in the ab¬ 
sence of relaxation. We find that in this limit the con¬ 
tributions to W2,W3) involving nonresonant 

transitions scale as F^. For resonant transitions, there 
are two cases that require further attention: (i) Taking 
Wcomb to be a possible frequency combination appear¬ 
ing in the expression in Eqs. (irara . resonant transi¬ 
tions (real or virtual) occur as \hujcomh\ = 2 |/r|. Then the 
function 77 ^ or becomes 'H^iiiUcomh + oc F“^ 

or I^(wcomb + i^'^Fe) OC F“^ respectively, and then 
^(3),dabc ^ p^p-1. i^g limit depends on the sequence of 
limits of Fi —>■ 0 and Fe —>■ 0, and so there seems to 
be no single well-defined relaxation free limit within this 
phenomenological theory, (ii) For some Wcomb = 0 , there 
can be divergences that occur at real frequencies in the 
absence of relaxation; below we discuss the behavior of 
^(3);dahc these divergences by considering the fre¬ 
quencies in the neighborhood of some of them. 

B. Finite temperature 

In calculating the response of a system to optical radi¬ 
ation, two effects of the temperature are usually consid¬ 
ered: its role in establishing the initial electron distribu¬ 
tion, and how it affects relaxation rates. In this work, the 
latter is implicit in our choice of relaxation rates. In our 
perturbative calculation, the former can be taken into ac¬ 
count in the following simple way: Explicitly displaying 
the chemical potential and temperature dependence, we 
write nsfc(/r, r) for the electron distribution at equilib¬ 
rium, and {fj,, T) for the nonlinear conductivity. By 
using 

^00 Q 

nsk{^J-,T) = J dxFfj,{x,T)—nskix,0) ( 33 ) 

with F^{x,T) = [I -|- the conductivity at fi¬ 

nite temperature can be related to the zero temperature 
conductivity via 


^(3);dabc(^^,^2,W3) 


io’s^o 

6\^\h^UJiUJ2UJ3 


+ 0 {u;-^), ( 32 ) 


in which the leading term is exactly the same as the con¬ 
tribution proportional to the 5 i term, and which agrees 
with the relaxation free, semiclassical calculation^^ in¬ 
volving only intraband motion. While the physically ap¬ 
propriate result of purely intraband, semiclassical motion 
is recovered in this limit as it should be, the connection 
to formally intraband, interband, and mixed responses in 
tT(3);‘'“^^(a;i,a;2,a;3) is less than direct. 

Erom Eqs. dSOllSIl) we can also study more generally 
the limits as the relaxation rates are allowed to vanish. 
Here we discuss the simple case where the intraband and 
interband relaxation rates are the same for all orders, 
but perhaps different than each other: F) = F^ and 
Fe ^ = Fg. We find that as F^/g —>• 0 we recover from 


/ oo o 

dxF^{x,T)—a^^'> ix,0) 

/ OO 

dxF^ (x, T) [I - (cr, T)] a( 3 ) (x, 0 ). ( 34 ) 

-OO 

Here the second line is obtained by using the partial in¬ 
tegration and the condition (j^^\x —>■ ±oo, 0 ) = 0 . Be¬ 
cause F^{x,T)[l — Ff^{x,T)\ is a pulse function located at 
X = fi with a width of the order of the thermal energy, the 
conductivity at finite temperature T can be obtained by 
averaging the zero temperature values over the chemical 
potential in an energy window with a width of the order 
of magnitude of the thermal energy. In a case where the 
chemical potential /i and the frequencies {uji} are chosen 
to be away from resonant transitions, the conductivity 
is a smooth function around /r. Considering that the 
thermal energy ksT is only about ^ 25.8 meV at room 
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temperature, the conductivity at room temperature is 
close to the value at zero temperature away from reso¬ 
nant transitions. However, around resonant transitions 
where the conductivity diverges, the effects of finite tem¬ 
perature can be important. In Eqs. (fT^ to (1261) . the 
chemical potential appears in the functions Gfj., 'Hfi, and 
I. in (T^ and cts, and as |/r| ^ in (tq- Therefore, the con¬ 
ductivity at finite temperature is determined by applying 
Eq. (l34t to these quantities. The temperature effects on 
the contributions due to the functions and 

are discussed in Appendix [BJ 

Note that the treatment of the ac term requires partic¬ 
ular care, because at finite temperature there are always 
electrons initially near the Dirac points, and they will 
lead to the same prediction for divergent response that 
Eq. (1311) indicates for electrons near the Dirac points at 
zero temperature in an undoped sample. To show this 
explicitly, from Eq. (l34l) . |^|~ should be replaced by 

/ oo 1 

dxF^{x,T)[l-F^{x,T)]— (35) 

-oo HI 

However, this diverges due to the singularity of the in¬ 
tegrand at x = 0. Based on Eq. (EH where this term is 
nonzero only at T^'^^ ^ the divergence shows that 
either the perturbation theory or the assumption of un¬ 
equal intraband and interband relaxation times in un¬ 
doped graphene is not adequate, and more realistic treat¬ 
ments of the scattering and temperature are required. 
Nonetheless, from a full numerical solution of Eq. m and 
we find that contributions from the |/r| ^ term only 
give a small contribution to the total conductivity at fi¬ 
nite temperature. Thus, at least at the level of the full 
SBE, whatever the final description of relaxation yields 
for the term it will not lead to significant contri¬ 

butions. So for our finite temperature calculations we 
somewhat arbitrarily take 

A ^ ^ . (36) 

ImI ^/JF+Jk^ 


C. Pulse response 

Because most nonlinear experiments are carried out 
using laser pulses, the optical response close to the diver¬ 
gences mentioned above is determined by the pulse shape. 
Except for the divergences just discussed, the in- 

elusion of the relaxation parameters moves the di¬ 
vergent frequencies off the real axis. Yet it is necessary to 
investigate the pulse effects when the energy broadening 
of the pulse is larger than the broadening characterized 
by those relaxation parameters. For a field associated 
with pulses of a fixed polarization 

E{t) = Y,Eu:,Pc..{t)e~'"'* (37) 


with the time domain envelope function Puiiit), the 
Fourier transform is 

E{u}) = ^ E^.P^.(uj - uJi) (38) 

i 

with the frequency domain envelope function 

Puiiiuj) = J 


The third order current in both time and frequency do¬ 
main can be written as 


J^{t) = 

Imn 


Imn 


with 


f dSidS 

m dSr, 


(39) 


j (2^)3 ' 

XaF),d.o-bc^^^ _|_ SijUJm + + ^n) , 


and 


= J ^j^PMP^^{6m)Po.Ad -5i- S^) 

_|_ Si,UJm + + S — 6l — 6m) -(40) 

We will be particularly interested in two special cases: 
(i) For 6i sufficiently small and ci-(3)’da.bc sufficiently slowly 
varying in its frequency dependence so that 

(41) 

over the frequency components of the envelope functions, 
the current response is given by 

(42) 

For a Gaussian pulse Pui^it) = e“‘ which gives 
P^,{uj) = we get 

^ , (43) 

.(44) 

with A“^ = A)"^-|-A“^-|-A“^. In this case, the generated 
currents are also Gaussian in their time and frequency 
dependence. 

(ii) For singular behavior 


0'^^'^’‘^°'^c(^UJl+Sl,UJm+6m, bOn+6n) 


ir^i3ydabc^u;i,u;m,u;n) 

6l + 6m + + *7 

(45) 







where 7 contains contributions from the relaxation pa¬ 
rameters, the optical coefficient (t) satisfies 

The solution of this equation is 

= f dre^^ 

J —00 

xpcji (t + t)p^^ {t + t)p^^ {t + t) ,(47) 

and 



ndabc 


id) = 




S + ij 

X J dte^‘^'^puii it)puj^ it)Puj„ (t) . (48) 


For a Gaussian pulse, we get 

^dabc (^\ _ (3);c?a6c/ /(t/A, A 7 ) 

= V [UJl,UJm,UJn) - 






S + A 


where I{x,y) = ^ye [1-I-Erf (a; — y/ 2 )], and 

Erf(a;) is the error function. In the absence of relaxation, 
we have lim _ .0fA/2[l -b Erf(t/A)], which is 

a constant as t —>■ 00 . This means that the current 

is nonzero even after the optical pulses have passed, in¬ 
dicating that current injection has occurred. Eor finite 
7 i (t) at t —)> 00 is zero, but the injected cur¬ 

rent can still persist for some time. Fig. [5] shows the 
dependence of the current response on the pulse width. 
For a very long pulse, 7 A ^ 1, the current response has 
a shape that is nearly Gaussian; however, for 7 A < 1, 
when the energy broadening of the pulse is larger than 
the relaxation rate, the current response obviously de¬ 
viates from Gaussian shape, and can last long after the 
excitation pulses are passed. 


III. THIRD ORDER OPTICAL 
NONLINEARITIES 

To illustrate how relaxation affects the third order op¬ 
tical nonlinearities, in the sample calculations presented 
below we assume equal relaxation rates for all orders of 
response, putting F^”^ = F^ and ri”^ = Fg, and consider 
four sets of parameters: (a) F^ = Fg = 0 , (b) Fj = Fg = 
33 meV, (c) F^ = 65 meV and Fg = 0.5 meV, which 
are parameters used by Gu et alA (d) F^ = 0.5 meV and 
Fg = 65 meV. We define set (a) by the limit F^ = Fg —>• 0, 
which recovers our relaxation free calculationi^ 


FIG. 2. (color online) Time evolution of I{t/A, 7 A) for 7 A = 
0.2 (red), 1 (blue), and 5 (green). The Gaussian pulse is 
plotted as black curve. 


A. Third harmonic generation 

For monochromatic incident light with frequency w, 
light is nonlinearly generated to lowest order at the 
third harmonic frequency 3a; and at the fundamental 
frequency ui. The first is described by the conductiv¬ 
ity the second corresponds to Kerr 

effects and two photon absorption, both described by 
CT(3);‘^“^°(_a’, O’, O’), and can be considered as a nonlinear 
correction to the linear optical response. In this section, 
we consider THG. 

The conductivity tensor for THG only has one inde¬ 
pendent component 

A^)’^^yy{uj,uj,uj) = A^^’^y^yiuj,uj,uj) 

= A^^’^yy^{uj,uj,uj) = (T(3);^^“^(u;,a,w)/3. (49) 

The induced current responsible for the THG is 
4 hg(^) =e"'^‘"V(3)™(cu,a,w)F;<,F;,,-F;,,+c.c. (50) 

In Fig. ([3]) we give the result for at 

\p\ = 0.3 eV for zero and room temperature. We first 
look at the results for zero temperature. The relaxation- 
free results are given as the thick (red) curves in Fig. [3] 
(a). This figure shows the step function of the real parts 
and the logarithmic divergence of the imaginary parts 
at three resonant photon energies hu> = 0.2, 0.3, and 
0.6 eV, which correspond to the resonant transitions for 
which the chemical potential gap 2\p\ matches the en¬ 
ergies of three photons, two photons, and one photon, 
respectively!^ With relaxation included, the conductiv¬ 
ity is a smooth function of a>, and plotted in Fig. [3] (b), 
(c), and (d). Some common effects induced by the relax¬ 
ations are shown: (i) the divergent peaks of the imagi¬ 
nary parts of the conductivity in Fig. [3] (a) become finite 
and broadened, (ii) the step functions of the real parts 
become continuous, (iii) the real parts become finite as 
huj < 2|/i|/3 = 0.2 eV, and increase rapidly with de¬ 
creasing frequency. They receive contributions not only 
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huj (eV) huj (eV) 

FIG. 3. (color online) Spectra of at zero (thick red curves) and room (thin bine curves) temperatures for 

different relaxation parameters: (a) Fi = Fe = 0, (b) Fi = Fe = 33 meV, (c) Fi = 65 meV and Fe = 0.5 meV, (d) Fi = 0.5 meV 
and Fe = 65 meV. The real (imaginary) parts of the conductivity are given by the solid (dashed) curves; we have taken 
1^1 = 0.3 eV. The insets focus on resnlts in the region [0.5, 1] eV. In (c) the fine strncture in region [0.198,0.202] eV is also 
displayed. 


from intraband transitions, describing Drude-like effects, 
but also from the interband transitions due to the lin¬ 
ear dispersion relation of graphene (for example, see the 
prefactor in Eq. (15^ ). 

To illustrate the dominant features in these fine struc¬ 
tures, we can analytically expand the coefficients of the 
functions and the term in the conduc¬ 

tivity, for small relaxation parameters ri^e/(?iw) <C 1, to 
write 


{3)-,xxyy 


(u, w, w) = -b tTg 


(3);xxyy 


(w) 


(3)\xxyy 


(w), 


(51) 


with 


(3)-,xxyy 


(w) 


ia3 


ug^j^ihu -b iVe) 


U4(hujf 

-640^(2^0; -b iTe) + 450^(3^0; -b zFe)] , (52) 


{3);xxyy 


Uj) 


CT3 


— 8'Hy,{2huj -b iTe) 


h 36(Sw)'^ - 
-bl77^^(3/kn -b iTe) T 3wX^(3Ikn -b iTe) 


^^3);x.yy(^)^_(r_r^)2 2^a3 


27{Huj)^\fj,\ 


In the relaxation-free limit as Te,i —>■ 0, —>■ 0 

and recovers the results of our previous worki^ 

However, the relaxation-free limit of strongly 

depends on the details of the chemical potential and re¬ 
laxation parameters; this is the contribution to 
from the general term discussed earlier in Eq. (ED, which 
is problematic unless E^ = Eg. Eor doped graphene where 
/i is finite, goes to zero with decreasing relax¬ 

ation parameters; for graphene that is undoped or at low 
doping, a more sophisticated treatment is in order, as 
discussed in Section III Al Eor the limit Ej_e,|/.t| "C fiw, 
the THG coefficient is approximated as 


^(Shxxyyi'^^ 


-ia3 \ 16(E,-Ee)" 

72{hjjj)'^ ^ 3fiw|^| 


(53) 


The term proportional to |^| ^ did not arise in our pre¬ 
vious calculation®, where we assumed that E^^e —> 0 
faster than /r —>■ 0. Deferring the treatment of small 
doping to later studies, we focus here on graphene with 
large enough chemical potential that (uj) does 

not make a significant contribution to the full third har¬ 
monic conductivity. 
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At room temperature, the conductivities for different 
relaxation parameters look very similar to each other, 
and the fine structures caused by the resonant transitions 
are smeared out. This can be understood by the results 
in Appendix |B] temperature affects the conductivity by 
smearing and lowering the peaks caused by functions C/, 
H, and I, which has an effect similar to increasing the 
value of Tg. If we increase each Te by the thermal energy 
of room temperature, the values of these new Tg in the 
four cases presented in Fig. [3] are close, and it is not 
surprising that we get similar room temperature results. 


B. Kerr effects and two photon absorption 


We now turn to the light nonlinearity generated at the 
same frequency w of the incident light. Taking £1^ = 


xE^ + yE ^, we write 1 consider 


Eui — Ei^ 


COS(j) 

sin 


,e-u, = e: 


COS(j) 

sin 


• (54) 


The nonlinear response at frequency ui is then given by 


J(3)(cu) = 3E^\E^\^ (a;,c., -a;) 

UJ, -Uj)i sin(2(/.) sin 9 


— cose 


• (55) 


For linearly polarized light (0 = 0), the second term van¬ 
ishes; the current from the first term has the same po¬ 
larization as the incident field, and gives an intensity de¬ 
pendent correction of the linear conductivity a^g{uj) = 
+ with 


aM = 3ct(3)™(c^,cc,-c.)|£;^|2 . (56) 

An effective nonlinear susceptibility can be 
introduced^^— Xni{^) = o-„;(a;)/(—fweodgr), where 
the effective thickness of graphene single layer dgr is 
taken to be S.SA^ from this an effective nonlinear 
refractive index n 2 and nonlinear loss /3 tpa can be 
extracted. In general, has both real part and 

imaginary parts, and the calculation of n2 and /3 tpa 
should follow the results of del Corso and Soles 

In the limit of no relaxation, we showed earlier— that 
cr(3);‘^“^°(_a;, cu, w) has many divergences, and its behav¬ 
ior in the neighborhood of equal frequencies can be writ¬ 
ten as 

a(3);''“''"(-w,a; + 5i,w + d2) 

(w; <52) 61 ) 

S 162 ^ h ^ I 2 

+r3"“'^(c.;<5i,52), (57) 

where and are all smooth functions 

of 5i and <52. The strength of the singularities is deter¬ 
mined by and which are real functions and 


are only nonzero when the photon energy is greater than 
that for the onset of one-photon absorption {Hlj > 2|p,|). 
For fixed photon energy fku, the appearance of the di¬ 
vergence as y decreases from 2 \y\ > Iia;to2|^| < huj 
indicates that it is associated with the existence of elec¬ 
trons (holes) at the k where one-photon absorption is 
possible. Physically, at these k the third order correction 
to one-photon absorption would lead to the perturbative 
description of the saturation, but in the absence of relax¬ 
ation that correction diverges, as it would for an inhomo- 
geneously broadened collection of two-level systems. At 
zero temperature, the sharp Fermi surface can strictly 
exclude such electrons (holes) for 2\y\ > huj. However, at 
finite temperature thermal fluctuations will always place 
some electrons (holes) at k where one-photon absorption 
can occur, and so the divergence in the third order re¬ 
sponse will exist for any photon energy. 

Including relaxation parameters F) = F^ and Fg = 
Fg, -I- 5i, w -I- ^2) includes terms that are 

proportional to ti = (5i -I- iFi)”^, t 2 = {S 2 + 
ta = (A -I- <52 + 2iTe)~^, tita, and t 2 ta- As any of these 
three quantities 5i -I- iFi, <52 -f fF^, or <5i -I- ^2 -f 2irg goes 
to zero, diverges. But for nonzero F^ 

and Fg, + + is a smooth function 

of real uj, 5i, and ^2. For a pulse response when the 
energy broadening of the pulse is less than the relaxation 
energies, it is reasonable to set 5i = ^2 = 0, and then the 
conductivity can be written as 


^{S)\dahc 


{-UJ,UJ,UJ) 


£dahc(^) ^ 

F.Fg f) 

/:f^^(cc;F„Fg) 

Fg 

+/:f''^(a;;F„Fg), (58) 


with 


"^(c^) = las 


12(fia;)' 


[g^i^+) + g ^{^.)], (59) 


/:f''=(cc;Fg) = 




12{huj)^ iz??. 
5i9l - 302 


2dl 




2hjjj 
^ L 




4r2 


(60) 


£r^(b.;r„rg) = 


rjg 




(s{hujYv+ 


\y\{Ay? + Tl)v+ 


|f(Ao - Ai)g^{iTe) 

-|- Aq H—I ,( 6I) 


The full expression of £4“^^(a;; F^, Fg) is complicated; we 
can achieve a good approximation by setting F^ = 0, for 
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TvjJ (eV) 



hu (eV) 



hw (eV) 



hjj (eV) 


FIG. 4. (color online) Spectra of at zero (thick red curves) and room (thin blue curves) temperatures with 

different relaxation parameters: (a) Fi = Fe = 0, (b) Fi = Fe = 33 meV, (c) Fi = 65 meV and Fe = 0.5 meV, (d) Fi = 0.5 meV 
and Fe = 65 meV. The chemical potential is |/r| = 0.3 eV. The real (imaginary) parts of the conductivity are given by the solid 
(dashed) curves. The result in region [0.1, 0.55] eV of (c) is highlighted in the inset. In figure (a) the real part (solid curve) of 
oj, cu) diverges for all fio; > 2 |p|. 


which 

£f'’=(w;0,re) 

= - ^o)eM(*re) + 4(3Ai + A^)g^{d+) 

— (4^1 + 5Ao)l/^(i9_) — 8(Ai + Aq)Q fj^(2huj + iVg) 

^2 _ 4 2-1 

(Ai + AAo)^+ + I. ( 62 ) 

In these expressions, we used = i/ku + iVe and iy± = 

±haj + iTi. 

In Fig. m the photon energy dependence of 
tT(3);^“^‘=(_a;,w, w) is plotted for different relaxation 
parameters, with chemical potential |/r| = 0.3 eV 

at zero and room temperatures. Fig. |4] (a) gives 
the relaxation-free calculation, which is done as 

lim Three regimes are apparent: 

r ie —^0 

( 1 ) huj < \fi\, in which both one- and two- photon ab¬ 
sorption are absent, and the real part of the conductiv¬ 
ity is zero. The imaginary part at low photon energy 
scales as {huj)~^. At huj = |/r|, the real part shows a step 
function, while the imaginary part shows a logarithmic 
divergence. (2) \fi\ < Hlo < 2|p,|, in which two-photon 
absorption is present but one-photon absorption is still 


absent. The real part of the conductivity here scales as 
Around huj = 2|/r|, the imaginary part shows 
a divergence {huj — 2|/r|)“^. For frequencies satisfying 
huj < 2\fi\, if the graphene is subject to a Gaussian pulse 
sufficiently narrow in frequency, the nonlinear current in¬ 
duced will still have a shape that is approximately Gaus¬ 
sian, and characterizing the nonlinear response to a pulse 
by Eq. 1(^)1 makes sense. (3) fiw > 2|/r|, where both two- 
and one-photon absorption are present. The imaginary 
part of the conductivity diverges as (Hlo — around 

hiu = 21 ^ 1 , and the real part diverges for the entire region 
tujj > 2|/i|. At finite temperature and in the absence of 
relaxation, the real part diverges for any photon energy 
huj. As we discussed after Eq. m, the divergence of the 
real part of oj, oj) is induced by the existence 

of electrons (holes) at the k where one-photon absorption 
occurs; at zero temperature, these electrons (holes) only 
exist when the chemical potential |/i| < huj/2, while at 
finite temperature, they exist at any chemical potential 
due to thermal fluctuations. 

In Fig. 0] (b), (c), and (d), we present the results for 
the same relaxation parameters as those adopted in the 
THG calculation. Relaxation affects the conductivity in 
a complex way, but there are some qualitative features 
that can be identified: (i) In the neighborhood of the 


4fj,{hu}y 

( r ?2 - 4 ^ 2)2 
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divergences that arise in the relaxation-free calculation, 
including the divergent regime ftu! > 2 \fx\ and the special 
frequency hu = |/r|, both the real and imaginary parts of 
the conductivity are lowered and are everywhere finite. 
In Fig.|4](b), (c), and (d), we find that a larger Fe gives 
lower and broader peaks at huj = \fi\ and 2|/r|. (ii) For the 
relaxation parameters used here, the real part of the con¬ 
ductivity is negative for huj > 2\fi\. Because of the pres¬ 
ence in this frequency range of one-photon absorption, 
which is always positive, the two-photon absorption pro¬ 
cesses indicated by the real part of w, w, uj) can 

be understood as a correction to the simple linear pre¬ 
diction of the absorption. In fact, we can find a range of 
electric fields large enough so that cr^g(u>) is negative; for 
a field anywhere near or above this strength the perturba¬ 
tive result is naturally suspect, (hi) Even for frequencies 
in the range |/r| < hoj < 2 |/r|, where only two-photon ab¬ 
sorption is present in the absence of relaxation, the real 
part of the nonlinear conductivity can 

be negative. Yet in the presence of relaxation the lin¬ 
ear conductivity acquires a real part in this 

frequency range, and the real part of (Tgg(w), for exam¬ 
ple, is always positive for small enough electric field am¬ 
plitudes, indicating absorption. However, these results 
indicate the sensitivity to the relaxation parameters of 
both the third order conductivity, and its interplay with 
the first order conductivity, and a more sophisticated de¬ 
scription of the scattering is clearly in order. 

At room temperature, the peaks or divergences are fur¬ 
ther broadened. For a given frequency w, the regime 
huj > 2\fj,\ always contributes to a finite temperature cal¬ 
culation due to the average over the chemical potential. 
The absolute values of the real part of the conductivity 
in the regime fko < 2 |/x| also significantly increase. 


C. Two-color coherent current injection 

Now we turn to two-color coherent current injec¬ 
tion, with frequencies wi = W 2 = —w and wa = 
2a;. In the relaxation-free calculation, the conductivity 
CT(3);‘^“^°(_a’,—O', 2a;) diverges, and it is the divergence 
that describes the current injection. In fact, in the neigh¬ 
borhood of these frequencies, the conductivity can be 
written as 

^(3);dabc(_^^ -w, 2w + 5 lu ) = , 

60UJ 

(63) 

where the injected current is determined by a well- 
behaved function (uj), and tT^^’‘^“^^(a;) is a smooth 

function of Sco. With the inclusion of relaxation, the 
conductivity a, —a, 2 w) itself is well behaved. 

The divergence term in the relaxation-free limit becomes 
a term similar to the right hand side of Eq. (HSl) . To check 
whether F^ and Fg have the same importance for the in¬ 
jection, we give the pulse calculations of 
(see Eq. 133) in Fig. [S]for different relaxation parameters. 
After the laser pulse, the current response persists for 



FIG. 5. (color online) Time evolution of Im[C'f'^^f.,^_ 2 i,;(f)] 
for Gaussian pulses with = 100 fs for different 

relaxation parameters (ri,re) = (0.2,10) meV (blue chain 
curve), (1.0, 10) meV (green dotted curve), and (10,10) meV 
(red dashed curve), (10, 1) meV, and (10,0.2) meV. The last 
three cases overlap with the red dashed curve. The pulse 
width corresponds to an energy broadening hA-i^ = 6.6 meV. 
Our calculations show Re[C'f.)ii(f.„^ 2 i.u(t)] is negligible on this 
scale. 


times associated with Fj, showing that the contribution 
from Fi dominates the relaxation of the injected current, 
as might be expected. To highlight this, we write 


a 


(3)‘,dabc 


(—UJ + (5i, —UJ + 62 , 2uj + 63 ) 


■^da 6 c(^) 


3(i5i S 2 S 3 ih ^Fi) 


{3);dabc 


(uj) , (64) 


From Eq. m, only terms including (5i -|- <52 + <^3 + 
iFi)”^ contribute to (a). By writing crf“j'‘^(a) = 
Hr]^abcI qY^ the first term in Eq. (1641) becomes 


_leV_ 

3(Mi -I- hS 2 + M 3 -b iFi) 


dabc{^ 
^inj 


'io'Z:T[uj) 


(65) 


with 


a. 


ia 3 J 2A3-I-A0 

, ^3 + ^0 


(M)3 

6reA3 -|- 5FiAo 


Im[£/^(2/ia; -b iFe)] 

^ Re[gf,(huj + iTe)] 


4:(Huj) 

3re(Ao -b A 3 ) -b FiAg 
4(M;)4 
ri(A 3 -b Aq) 


Re[gf^( 2 huj -b iFe)] 


4(M;)3 

ri(A3 -b 2Ao) 

4:(tlUj) 


'Hfj.(—2hjj + iTe) 

3 'H^(/ia;-b fFe) I , ( 66 ) 
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FIG. 6. (color online) Spectra of (a, b) and (c, d) for different relaxation parameters Vi and Fe at temperature 

T = 0: Black solid curves; Fi = Fe = 0; Red dotted curves: Fi = Fe = 33 meV; green chain curves: Fi = 65 meV, Fe = 0.5 meV; 
Blue dashed curves; Fi = 0.5 meV, Fe = 65 meV. The chemical potential is |/i| = 0.3 eV. The result of (to)] at 300 K 

is shown in the inset of (b). 


where terms proportional to are neglected. As 

FijFe —>• 0, only terms involving Im[t/^(a;)] oc 9{\uj\ — 2\fi\) 
remain; is a pure imaginary quantity, and is 

consistent with our previous worki^ With relaxation in¬ 
cluded, terms involving Re[C/^] and T-L^ appear. 

In Fig. [6] (a) and (b) we plot the photon energy de¬ 
pendence of for different relaxation parameters 

at zero and room temperature. The real part is much 
smaller than the imaginary part. At zero temperature, 
we see in Fig. [5] (b) that for Fg = 0.5 meV there are 
fine structures in the spectrum of Im[crf^j““(a;)] around 
fiw = \^\ and hw = 2|/x|; it is due to the T-L^ terms. As 
can be seen in the inset of Fig. [5] (b), finite temperature 
and finite Fg lead to similar broadening and lowering of 
the peaks. 

We also plot in Fig.|n](c) and (d). The ampli¬ 

tude of is of the same order of magnitude as that of 
However, due to the prefactor leV/(3Fi), which 
relates to its contribution to w, —w, 2a;) 

in Eq. (1^ . dominates for small F^, usually taken 

to be a few tens of meV. 


D. Parametric frequency conversion 

Third order nonlinearities can lead to the appearance 
of new frequencies via parametric frequency conversion, 
which is described by a;s,a;p,a;^). Here ujp is 

the frequency of a strong pump field, and ujg is the signal 
frequency converted by interaction with the pump to an 
idler frequency iOi = 2ujp — ujs . Possible resonant transi¬ 
tions occur as any of the frequencies |a;p|, |a;s|, jojp — oJsI, 
|2a;p|, or |2a;p — a;s| equal 2\^\/h. In Fig. 0 we plot the 
dependence of on uJb for different 

relaxation parameters at hwp = 0.8 eV and |/i| = 0.5 eV. 

At zero temperature, the calculations show peaks/step 
functions for resonant transitions at Hujsi = 2 Hujp—2\fj,\ = 
0.6 eV, fujJs 2 = 2|/r| = 1.0 eV, or hwss = 2|^| -|- fiwp = 
1.8 eV, both with and without the inclusion of relaxation. 
Around these resonant transitions, the behavior of the 
conductivity can be analyzed as following: 

( 1 ) Around ujg = uisi + Sui, the idler photon energy 
huii = 21^1 — hSu! is close to the onset of the one-photon 
absorption. By taking -d = huji + iTe = 2|/r| — hSuj + iTg, 
the conductivity as (5a; —>■ 0 is determined by functions 
t/p(i?) ln(/i(5a; - IFg), {HSoj - ire)“iFi, and 
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FIG. 7. (color online) Spectra of {—LUa,ijjp,u!p) for different relaxation parameters at zero (thick red curves) and room 

temperature (thin blue curves). Solid (dashed) curves give their real (imaginary) parts separately. In the calculation, the 
chemical potential is |/i| = 0.5 eV, and the pump photon energy is hup — 0.8 eV. 


2)i(i?) {hSu! — iTe)~^Ti. In the relaxation-free limit, 
only Gfj.i'd) contributes a logarithmic divergence to the 
imaginary part, and a step change in the real part [in 
Fig. 13 (a)] for nonzero Soj. With the inclusion of relax¬ 
ation, we can distinguish three different types of qualita¬ 
tive behavior, shown in Fig. |3(b) - (d), based on the rel¬ 
ative magnitude of F^ and Fg: (b) Fj = Fg, all functions 
contribute; (c) Fg <C Fj, I^(il) dominates; (d) Fg » F^, 
where for the values chosen the relaxation is large enough 
to smear out these resonances. 

(2) Around Ws = 0 Js 2 + Suj, the signal frequency is 
close to the onset of the one-photon absorption. For 
non-resonant transitions in a usual semiconductor, 
and Wsi are interchangeable frequencies to give the same 
conductivity of parametric frequency conversion)^ here 
in graphene they yield asymmetric peaks because the res¬ 
onant transitions dominate. In the limit of no relaxation, 
the conductivity shows a logarithmic divergence that is 
easily smeared out by the inclusion of small relaxation 
parameters. 

(3) Around Wg = Wgs-l-i^a;: By taking d = —hu!s + hujp+ 
jFg = — 2|/r| — h6uj + zFg, the conductivity as (5w —>■ 0 
is determined by functions Qp,('d) ^ ln(Ma; — fFg) and 
'H^('d) ^ (Mw — frg)“^ri. In the limit of no relaxation, 
Gp. gives a small peak. For Fg <C Fj, the peak from 'Hp 
is stronger but very narrow. 


(4) At finite temperature there is a further smearing 
of the peaks around the resonances, as described in Ap¬ 
pendix |b] 

Besides these resonant transitions, two singularities are 
apparent: (i) the singularity around Wg = 2a;p, which 
corresponds to two-color coherent current injection. The 
singularity is not determined by the behavior of the func¬ 
tions Gi 'H, T, but by the coefficients that premultiply 
them in Eqs. (IT^ to (l 26 ll . Around Wg = 2 a;p we put 
Wg = 2u}p + Slu and find 


a(3);""""(-2u;p + dw,wp,wp) 


leV 




P 

{3)\xxxx 


hSuj + iVi 


/ \ , (3)\xxxx / \ 

{u}p) + ay [ujp), 


an equation similar to Eq. (IMl) discussed in Section ITlI Cl 
Two-color coherent current injection requires both one- 
photon absorption (for Wg) and two photon absorption 
(forwp), ie., hujp > \^\. The parameters we have adopted 
in Eig. 13 fulfill this criterion, and thus the singularity ap¬ 
pears. Finite temperatures do not qualitatively affect 
this singularity because it is not related to the chemical 
potential, (ii) The strong response around Wg = uip is 
related to the third order correction to one-photon ab¬ 
sorption, and only appears at finite temperature. Since 
1^1 < hujp < 2|/r|, at zero temperature only two-photon 
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absorption is present and there is no one-photon absorp¬ 
tion. However, at finite temperature thermal fluctuations 
will place electrons where one-photon absorption can oc¬ 
cur, and the third order correction to that will lead, in 
the absence of relaxation, to a divergent result as dis¬ 
cussed in Section IIIIBl in the presence of relaxation the 
result will not be divergent but very large, describing the 
saturation of the one-photon absorption at the level of 
the third order response. 


E. Comparison between calculations and 
experiments 


Experiments have already extracted values of the ef¬ 
fective third order susceptibilities of two- 

photon absorption ;^^!^^ Kerr effects^ '^^i^^ and parame¬ 
ter frequency conversion^! at some photon energies. The 
nonlinear conductivities we have calculated here are re¬ 
lated to the effective susceptibility by^i^ 


(3):dabc/ \ 

Xj (wi,a;2,W3) 


— i{uJi L1J2 + W3)eo'^gr 


(67) 


We Hrst look at the THG, for which the experimental 
technique is perhaps the most mature, and the extracted 
values can likely be considered more reliable than those 
from other effects. For a reasonable chemical potential 
estimated from the sample preparation, the calculations 
without relaxation parameters^ yield theoretical results 
for the nonlinear conductivity about two orders of magni¬ 
tude smaller than the value extracted from experiments. 
Here we have found that calculations at finite tempera¬ 
ture for different sets of relaxation parameters (see the 
insets of Fig. El) are almost the same as calculations at 
zero temperature and neglecting relaxation. 

For Kerr effects, because of the existence of divergent 
terms in the expressions and the probably very low chem¬ 
ical potential in experiments, it is not surprising that we 
could fit the nonlinear susceptibility at one photon en¬ 
ergy by tuning the relaxation parameters. The compli¬ 
cated dependence is shown in Fig. E) As hu > 2|/r|, the 
nonlinear conductivity at both zero and room tempera¬ 
tures can vary many orders of magnitude, depending on 
the relaxation parameters adopted. 

For parametric frequency conversion observed in the 
experiment by Hendry et al.^ with parameters fiuip = 
1.31 eV, hws = 1-05 eV, and assuming a low chemical 
potential |/r| = 0.1 eV, we checked the dependence of 
the conductivity on the relaxation parameters T i and Tg 
in the range of [0,60] meV. We find the dependence is 
weak and the calculated values are still smaller than their 
claimed values by two orders of magnitude 4! 

Admittedly the measured effective susceptibilities for 
parametric frequency conversion, Kerr effects and two 
photon absorption, and THG show a strong dependence 
on the measurement method, light frequency, pulse du¬ 
ration, and perhaps sample preparation. Yet even taking 


this into account, the conclusion that the theoretical re¬ 
sults are about two orders of magnitude smaller than 
the measured results is inescapable. These discrepan¬ 
cies could arise for a number of reasons, including: (1) 
The samples in many experiments are not suspended 
graphene, but graphene on a substrate or in solution. 
Thus there may have been contributions to the optical 
nonlinearity from the interaction between the graphene 
sheet and its environment, which may be crucial consid¬ 
ering that graphene is a one-atom thick material. (2) 
Thermal effects^ caused by a high repetition rate of 
laser pulses, as used in Z-scan experiments, may play 
an important rolei^2i£i (3) Because of the zero gap of 
graphene and the intense laser beams used in exper¬ 
iments, saturation!! induced by one and/or two pho¬ 
ton absorption can make necessary a treatment more so¬ 
phisticated than that of perturbation theory. Zhang et 
al^ used the density matrix method to study four wave 
mixing in undoped graphene in the saturation regime, 
and found an effective xig about 10“^^ m^/V^, and 
decreasing with increasing light intensity. Additional 
calculations for different third-order nonlinear effects in 
graphene in the saturation regime are needed to assess 
the impact of saturation on the theoretical nonlineari¬ 
ties. (4) The calculation at the independent particle level, 
which works well as a starting point for most gapped 
semiconductors, may fail in graphene, and it may be nec¬ 
essary to do a more realistic calculation, including the full 
band structure, and the detailed effects of scattering and 
the electron-electron interactions. 


IV. DC CURRENT INDUCED SECOND ORDER 
NONLINEARITY 


We now turn to the limiting case where one of the 
electric fields is a dc field, taking a ;3 = 0. The calcula¬ 
tion of a; 2 , 0) from Eq. (fT^ includes a term 

proportional to 


1 


fiw3 + irf ’ ir 


-3(3) 


( 68 ) 


(3) 

Therefore a nonzero relaxation T) for the dc field is 
necessary to set up a steady state with a dc charge 
current in graphene. For other transitions included in 
(T(3);‘'“‘='’(a;i, 0 , 0 ^ 2 ) and W 2 ), it is not nec¬ 

essary to include relaxation associated with the dc field, 
because the dc field acts on the optical excitation with 
frequency 0 J 2 and wi -I- UJ 2 respectively; these only survive 
during the optical pulse. We list the relaxation parame¬ 
ters used in calculating the unsymmetrized conductivities 
in Table [H The third order conductivity of interest here, 
which we can refer to as the dc-induced second order op¬ 
tical conductivity, can be written as 


^(3);da.c(^^^^2,0) = ^ur'=(wi,a;2) 


dabc 
■<^E 


(wi, W 2 ) ■ 
(69) 
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unsymmetrized a 

relaxation parameters 

unsymmetrized a 

relaxation parameters 

ff(3);'^“'’"(cai,ca2,0) 

ff(®);'^'’“"(a;2,a;i,0) 

p(3) pdc 

p(3) pdc 

-L e — J- e 
p(l)2) _ pop 

_ pop 

ff(3)’‘^“"‘'(a;i,0,a;2) 

ff(3)'‘^'’“(a;i,0,a;2) 

ff(3)''^“'’(0,a;i,aj2) 

ff(3)''^"'’“(0,a;2,a;i) 

p(f?) _ pop 

rO) = pop 


TABLE I. Relaxation parameters used in the different processes associated with the dc current induced second order nonlinearity. 


The first term includes all contributions that di¬ 
verge as 1 /r^'^, which are only involved in calculating 
~( 3 );dabc(-^^^^ 2 , 0 ) and 5 ( 3 );dbac(-^ 2 ,a;i, 0 ); they both oc¬ 
cur with the dc charge current. Thus we can associate 
it with the dc current-induced second order conductivity, 
in that it is second order in the optical fields at wi and 
UJ 2 - The second term includes all other contributions, 
and we can associate it with a dc field-induced second 
order conductivity, which exists even for a gapped semi¬ 


conductor without doping. Examining Eq. (Ha , we see 
that W 2 ) is independent of T^'^ and T^'^, and it 

can be written as 




with 




Ai 

{huJ2) 


A^i, 


• + ■ 


A 9 . 


( 1 

1 ^1 

\hiui - 1 - tii02 + *r°p 

hU2 J 

( 1 

1 ^ 

V Swi - 1 - huj2 + iT°^ 

hioij 


Ai A 2 Ai A 2 + A^, f 1 1 \ 

{huj2Y fux}i + fux}2 + iT°fP \hwiiT°'^ tujj2 + iT°'^) 

As /I 1 M., . 


T- t ^ .pop ( 7-f T—) -I-/ia;2 + *r°P) 

fROi + riuj2 + li e \ nuji ruv2 J _ 

^tujj2 + iy°^ ^ -I-^ 


1 1 
tujji + hjjj2 + huji + fwj2 + ^ 

1 1 ' ’ ' 


x^(^huji + huj2 H“ 

° ( hwi + fr°p tujj2 + fr°p 


+ ^3 




(70) 


As we show below, the values of aj and ge are typically 
of the same order of magnitude; hence it is the value 
of TA that determines whether the dc-current induced 
second order conductivity or the dc-field induced second 
order conductivity makes the larger contribution to the 
dc-induced second order conductivity W 2 , 0 ). 

We can get a rough estimation of T^*'^ from the graphene 
mobility /Xm- The dc limit of the optical conductivity can 
be obtained from Eq. (HU as 




dgplMl 


(71) 


Hence, we get 


pdc ^ ^^F|e| 


(72) 


For a sample with mobility = 10^ cm^/(V-s) and 
chemical potential /i = 0.5 eV, T)*'^ is about 10 meV. We 
will see below that for samples with such mobilities the 
dc-current induced effects will typically dominate the dc 
induced second order conductivity. 


A. DC current induced SHG 


The connection between the mobility and conductivity 
can be written as (t^^^’^^(O) = Nf,\e\^rn with the carrier 

I 1^ 

density W = obtained from the linear dispersion. 


We first consider dc-induced second harmonic genera¬ 
tion, governed by w, 0). For monochromatic 

light at frequency w, the second order optically induced 
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FIG. 8. (color online) Spectra of (left column) and (right column) for different optical relaxation 

parameters at zero temperature (thick red curves) and at room temperatures (thin blue curves); |/r| = 0.5 eV. Solid (dashed) 
curves give their real ( imaginary) parts. In calculating ( 00 , 00 ), F))'^ = r°P and = 0. 


current is given by 

(73) 

There the two nonzero components are cr^HG(‘^) “ 
3cr(3);““i'i'(a;,w,0) and crsHG(w) = { 00 , 00 , 0 ). 

Correspondingly, each of them includes two parts: the 
dc-current induced second harmonic generation (CSHG) 
CTj ' and the dc-field induced second harmonic generation 


(EFISH) 

In Fig. [5] we plot the photon energy dependence of 
aj'^^£{oo,oo) for \fi\ = 0.5 eV and different values of opti¬ 
cal relaxation parameters r°P and r°^ at zero and room 
temperature. Two resonant peaks appear for both 
and one at hoo = \fj,\ and one at 2|/r|. The first 

corresponds to the second harmonic resonant with the 
onset of one-photon absorption, and the second to the 
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FIG. 9. (color online) Chemical potential |/r| dependence of aj°'’'^{—ujs,ujp) for different optical relaxation parameters at zero 
temperature (thick red curves) and at room temperatures (thin blue curves); hiOs = 0.58 eV and hwp = 1 eV. The real 
(imaginary) parts are given by solid (dashed) curves. The arrow in (a) indicates the vertical line there corresponding to the 
derivative of a <5 function. 


fundamental resonant with the onset of one-photon ab¬ 
sorption; the first peak leads to a higher response coeffi¬ 
cients than the second. In general, and are 

of the same order of magnitude. Therefore, in a high mo¬ 
bility graphene sample with a small Tf'^, the contribution 
of crj“^°(a;,a;) dominates because of the 

prefactor leV/T)*'^ (see Eq. (IMll ). 

From Eq. dZOl), we see that the first resonance in 
^<io&c(^) jg determined by 'H^(2huj + fr°P) and Ip,(2hw + 
jEgP); the other resonance is determined only by 
jEgP). Obviously, for both transitions smaller values of 
r°P result in a larger value and a sharper peak. At 
room temperature, these peaks are broadened and low¬ 
ered. The vertical line at /ku = \^\ in Fig. [5] (al) comes 
from lvci\H^{2huj + *r°P)] as r°P —>■ 0, which is propor¬ 
tional to 5{fu;j — 1^1); the peak in Fig. |5] (cl) shows the 
fine structure of Im|l^(2Scj -|- ir°P)] for r°P = 0.5 meV 
(see Appendix IbI). However, it is interesting to note that 
these two fine structures undergo important changes at 
room temperature: In Fig. [5] (al), we see that the first 
fine structure leads to a peak with broadened width; in 
Fig. |5] (cl), we see that the second fine structure leads 
to a sign change around tvjj = \^\ when the temperature 
increases from zero to room temperature. 


B. DC current induced difference frequency 

A counterpart of the third order parametric frequency 
conversion discussed in Section IIII PI is difference fre¬ 
quency generation which is, like second harmonic genera¬ 
tion, a second order nonlinear effect that can be induced 
in graphene when applying a dc field. With a strong 
pump at frequency Wp, difference frequency generation 
converts a signal frequency Ws to a new frequency ujp—ujs; 
the response is determined by Ws, Wp, 0). Sim¬ 

ilar to dc-induced second harmonic generation, there are 
current and electric field contributions to dc-induced dif¬ 
ference frequency generation. As we found in Section 
II V Al for dc-induced second harmonic generation, the cur¬ 
rent contribution should dominate the dc-induced differ¬ 
ence frequency generation in a high mobility sample. As 
an example, we plot the chemical potential dependence 
of Ws, Wp) for different optical relaxation pa¬ 

rameters in Fig. [^for Hujp = 1 eV (with a wavelength of 
about 1.24 /im) and hojs — 0.58 eV (with a wavelength 
of about 2.1 /im). For vanishing optical relaxation pa¬ 
rameters (r°P = r°P = 0), it is clear from Fig. [^(a) that 
there are three resonant transitions in the plotted chem¬ 
ical potential range: |/ti| = {hwp — huJs)/2 = 0.21 eV, 
|/t 2 | = hjJs/2 = 0.29 eV, and |/t 3 | = hjjp/2 = 0.5 eV. 
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Without relaxation, the imaginary part of the conduc¬ 
tivity is always zero except at these three resonant tran¬ 
sitions (shown as vertical lines): the first is given by 
lim - feug -I- id)] oc -£^5{hujp - hiOs - |mI): 

the other two are given by lim Swg -|- id)] oc 

5 —>^0 

6{hiUs — ImI) Etnd lim Im]Hp(?iaip -I- id)] oc 6{hiUp — jp,j). 
5 —>^0 

With finite relaxation rates or at finite temperature, 
the vertical lines are broadened to structures of finite 
strength and width. 

C. Comparison between calculations and 
experiments 


angle as Qcf) and 2>(j)] for the (001) face, the dependence 
is 0(j) and 4(/); while for the (110) face, the dependence 
becomes 0^, 2(j), and A(j). Therefore, from the azimuthal 
dependence of the SHG signal it might be possible to dis¬ 
tinguish the graphene CSHG (EFISH) from the interface 
contributions, for example by putting graphene on top of 
different Si02/Si structures, one with the (111) face of Si 
normal to the interface and one with the (001) face nor¬ 
mal. Because of the same origin of CSHG and EFISH in 
graphene, they would have the same angle dependence, 
so such experiments would not help to distinguish be¬ 
tween these different contributions from graphene; but 
for a heavily doped and high mobility graphene sample 
our calculations show that the CSHG should dominate. 


Bykov et ah^ observed that SHG radiation from a 
graphene/SiO2/Si(001) substrate strongly depends on 
the applied current density in the graphene layer, which 
is attributed to the CSHG effect of graphene. A sim¬ 
ilar structure was also studied by An et al.,— who 
could measure the radiation from different locations on 
the graphene sheet; they interpreted the result as EFISH, 
where the electric field is induced by current-associated 
trapped charge at the graphene/Si02 interface. Because 
of the interface contribution to the SHG radiation,— “— 
the contribution of the current related SHG from the 
graphene is hard to extract. 

The best way of measuring the dc-induced second or¬ 
der nonlinearity of graphene, without any background 
contribution from interface effects, would be to mount 
graphene in a symmetric structure; this can be difficult. 
However, within the framework of the experiments of 
the type that have already been done, we can suggest 
a strategy that might help identify the in-plane graphene 
GSHG(EFISH) by the azimuthal angle dependence of 
the generated signal. For linearly polarized light with 

and Fide = Adc Eq. ([TSD 

ysm0y \svcL(pj J 

becomes 


*^ shg (^) “ E^Edc^ai{uj) 


I cos(2(;i - c/)j) 
ysm{2(j) - (j)j) 


-|-]cri(a;) -I- 0-2 (w)] 



(74) 


We see that the Gartesian components of the induced 
current vary as cosinusoidal functions of 0(f> (that is, inde¬ 
pendent of (p) and 2(j). In a short-hand notation, we will 
characterize these as Qp and 2(j) dependences. In most 
experimentSfii^— the graphene sample is mounted on a 
Si 02 /Si substrate, where the interface between Si 02 and 
Si gives an interface-induced SHG and the bulk Si gives 
an electric quadrupole/magnetic dipole induced SHG. 
But for different crystal orientations of the Si substrate, 
the dependence of the combined interface and bulk con¬ 
tributions on azimuthal angle will be different^SS For the 
(111) face, the second harmonic radiation depends on the 


V. CONCLUSION AND DISCUSSION 

Perturbative analyses play a central role in nonlin¬ 
ear optics. Even when the electrons in a material are 
treated as independent, and relaxation is only described 
phenomenologically, the calculated response tensors that 
relate the induced polarization or current to powers of 
the applied fields indicate the nonlinear optical effects 
that are allowed, and point to where resonances can lead 
to an interesting dependence on time and frequency. Of¬ 
ten more sophisticated models of the electron dynamics 
are required, and sometimes the perturbative framework 
itself is insufficient to address the physics of interest. But 
even then these kinds of perturbative treatments provide 
a starting point for more realistic calculations. 

In this paper we have provided such a treatment of the 
nonlinear third order optical response of doped graphene, 
with the main goal of investigating the effects of phe¬ 
nomenological relaxation parameters, finite temperature, 
and laser pulse width on the induced currents. We fo¬ 
cused on the contributions of optical transitions around 
the Dirac points, where the widely used linear dispersion 
relation is a good approximation. By solving semiconduc¬ 
tor Bloch equations perturbatively, an analytic expres¬ 
sion for general third order conductivities was obtained 
at zero temperature, taking different relaxation parame¬ 
ters for interband and intraband optical transitions. The 
nonlinear conductivities at finite temperature were ob¬ 
tained by an appropriate integration over the chemical 
potential. The conductivities show a complicated depen¬ 
dence on photon energy, chemical potential, and the re¬ 
laxation parameters. 

Even with the inclusion of relaxation we found that the 
perturbative approach itself is problematic at vanishing 
chemical potential, as might be expected from a similar 
result in the semiclassical limit^^ except in the special 
case that either first- or third-order interband and intra¬ 
band relaxation rates are set equal. The perturbative ap¬ 
proach adopted is unproblematic for doped graphene at 
zero temperature, but is a concern at finite temperature, 
since thermal fluctuations always place some electrons or 
holes near the Dirac points. Yet numerical calculations of 
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the full semiconductor Bloch equations indicate that the 
contribution of such electrons to the full optical response 
is small, so this effect does not afflict our results. 

We discussed in detail different nonlinear effects, in¬ 
cluding third harmonic generation, Kerr effects and two- 
photon absorption, two-color coherent current injection, 
parametric frequency conversion, and dc-current and - 
field induced second harmonic generation and difference 
frequency generation. The interband relaxation gener¬ 
ally broadens and lowers the resonant peaks, while the 
intraband relaxation plays an important role in some of 
the effects, including two-color coherent current injection 
and the dc-current induced second order nonlinearities. 
At room temperature most of the resonant structures are 
smeared out. 

We also considered the response of graphene to laser 
pulses. The optical response depends in detail on the 
frequency width of the incident pulse and the frequency 
structure of the response tensors. The two natural limits 
are (1) when the frequency structure of the response co¬ 
efficients is rather flat on the scale of the frequency width 
of the incident pulse, and the induced current follows the 
injecting pulses, and (2) when there are divergences in 
the response coefficients at frequencies close to the real 
axis, as for two-color coherent control, the a dynamics 
associated directly with the relaxation processes. 

Comparison of our results with experiments is difficult, 
since in many of the reported experiments the graphene 
samples have not been characterized in the linear regime, 
and neither the relaxation parameters nor even the chem¬ 
ical potential have been identified. Results for some non¬ 
linear response coefficients, such as that describing the 
Kerr effect, are predicted by our calculations to be so 
sensitive to these parameters that we cannot hazard a 
comparison of theory to experiment. Yet our results for 
third harmonic generation and parametric frequency con¬ 
version are insensitive enough to these parameters that 
we can conclude our results are about two orders of mag¬ 
nitude smaller than those extracted from experiments^ 
even with the adoption of reasonable relaxation parame¬ 
ters. We speculated on the causes of this disagreement; it 
is of course early days for both detailed experimental and 
theoretical studies of such nonlinear effects in graphene. 
But these disagreements may persist, and a long and dif¬ 
ficult journey may be necessary to understand the details 
of the full nonlinear optical response of graphene. Even 
so we can expect that, as in the study of the nonlinear 
optical response of other materials, the kind of pertur¬ 
bative calculation we have presented here will provide a 
useful port of embarkation, paving the way to a better 
physical insight in the complex nonlinear optical response 
of graphene. 
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Appendix A: Perturbation solution of 
Semiconductor Bloch equation 


We describe the electronic states in graphene by the 
tight binding model employing carbon 2pz orbitals with 
only nearest neighbor coupling. Neglecting the overlap 
between different orbitals, the band structure that 
results is electron-hole symmetric with energies e+fc = 
—S-k where -l-(— ) is the band index for 7r*(7r) bands, 
and the Berry connections satisfy ^ssk — $.ssk = f'k and 
^ssk = ^ssk where s = ± is the band index and s is 
the index of the band that is not the s bandit Up to 
the dipole approximation of light-matter interaction, the 
SBE in Eq. ([T]) can be expanded as 


-dp. 


(") ( 4 .\ 
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dt 


- eE{t) ■ VkP. 


ssk 

(n-l) 


^ssk 
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= -^sekpi'iit) + ieE{t) ■ rkip^k “ ptk 
- eEit) ■ , (Al) 

with tk = £+fe — S-k- In Eq. (IA1|) the terms involving 
give the interband contribution, and the terms involving 
Vfc give the intraband contribution. Treating the elec¬ 
tric field term perturbatively, the first three terms are 
expanded as 


p^sLkit) = 

pf.l^ki*) = 
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with wq = W 2 + ^3 and lu = uii + uiq. By substituting 
the above expansion into Eq. we get the following 

equations for 

(1) The linear order terms are determined by 
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Here we have put Auk = n^k — n-k- The solutions are 
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which leads to the linear conductivity: 
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Here DsiS 2 fe the matrix elements of the velocity oper¬ 
ator, which satisfy v++fc = —v _^ and _^ = —v _|_fc 

in the tight binding model we have adopted. 

(2) The second order terms are determined by 

{■ffs - sek)V^f^'’‘'{uj2,uj3) = rl fiPi2"(‘^3) - v‘'^sk{u)3) 

(A6) 

urch 

The solutions are 
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Because the graphene crystal structure is centrosymmet- 
ric, its second order conductivity is zero. 

(3) The third order terms are determined by 
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The frequency dependence of { 101 ,^ 12 , 103 ) is im¬ 

plicit. The solutions can be written as 
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Here terms T’“^£ with i = I,-- - ,4 are related to the 
populations at band s and are given by 
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while the terms P“s£ with i = 5, • • • ,8 are related to the 
interband polarization and are given by 
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The unsymmetrized third order conductivity which fol¬ 
lows from these terms is 
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Then we find terms Si in Eq. (Ha are given by 
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for i = 1, 2,3,4 and 



for f = 5, 6 , 7, 8 . 

In this work, we only consider optical transitions 
around the Dirac points K = {bi + 2 b 2)/3 or K' = 
(62 + 26i)/3 with 61 and 62 the primitive reciprocal lat¬ 
tice vectors. These two Dirac cones are connected by 
the inversion symmetry, and they lead to the same con¬ 
tribution to the conductivities we consider, whether lin¬ 
ear or third order. In the following, we calculate the 
conductivity around K explicitly and get the total re¬ 
sults by considering both valley degeneracy = 2 and 
spin degeneracy gs = 2. Around the Dirac point K, 
we approximate each quantity up to its lowest order of 
fc — K: The electronic dispersion is £sK+k = shvpk, 
the velocity matrix elements are 'Uss(K'+fc) ~ svpk/k and 
'^ss{K+k) ~ isvpk X z/k, and the interband Berry con¬ 
nection is VK+k ~ fc X £/2fc2. The linear conductivity 
that results is given by Eq. (ED. 

In calculating the integrals over fc necessary to evaluate 
the third order conductivities, we use the relation 
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FIG. 10. (color online) t dependence of (a) + iF) and Gfj.-,T{£ + iF), (b) + iF) and + iT), (c) Xij,{e + iF) and 

Xfj,-T{e + iF). Solid (dashed) curves are for F = 0 (F = 0.03 eV), and thick red (thin blue) curves are for the function without 
(with) a temperature average. Other parameters are |/i| = 0.3 eV and T = 300 K. 


to expand the derivatives in 5^, and find that all the 
required integrations over k can be related to 


9sgv^ 

9s9v-^ 




1 , 
^0-3^0 , 


Using partial fractions and taking the integral over 0, we 
obtain all the Si given in Eq. (I191I26P 


Appendix B: Temperature Effects 


In this appendix we discuss how the temperature af¬ 
fects the contributions to the conductivities from (e + 
iTe), + iFe), and I^(e + zFe). 

(1) Gij,{i + i^e)' At finite temperature, this is replaced 
by + iFe) with 


Gfi-,T{e+i^e) = 13 / dxF^{x,T)[l-Ffj_{x,T)]Gx{e+iTe) ■ 

J —00 

(Bl) 

As Fe —t 0, t/^(e + iVf.) diverges logarithmically at e = 
±2|/r|, while Ge) is smooth. Both functions are 
smooth as e-|-fFe — 7 > 0. In Fig.[Tn](a) both Gfj.{i+i^e) and 
Gij,-Ti£ + *Fe) are plotted for Fe = 0 and Fg = 0.03 eV 
at |/r| = 0.3 eV and T = 300 K. Moving from zero to 
finite temperature, a finite T has an effect similar to the 
inclusion of relaxation: Both remove the singularity and 
broaden the peak and step function. 

( 2 ) + At finite temperature, this is replaced 


by + iFe) with 


H^;T(e+^Fe) = P / dxF^{x, T)[l-F^{x, T)]'Hx{e+iVe) ■ 

J —00 

(B2) 

At |/r| —>• zq with 2zo = e + *Fe for e > 0, -I- iFe) 

diverges as (|/x| — zq)”^- In the relaxation free limit for 
nonzero e we can write 


1 


|/r| - e - IF, 


re -*:0 p -- 1 _ ^ (-23) 


'ImI - ' 


where P. means the integration takes the principal part; 
thus the imaginary part of tends to a 5 function. 
However, both the real and imaginary parts of 'Hfj.-T are 
smooth for e > 0 or e = 0 and Fe 7 ^ 0. For small Fg, 
'Hii-T{i33e) oc InFe. 

In Fig. [in] (b) both + iTe) and 'H^;T(e + *Fe) are 
plotted for Fe = 10“^ eV and Fg = 0.03 eV a.t fj. = 0.3 eV 
and T — 300 K. The inclusion of finite temperature leads 
to a broadening of the (5-function-like imaginary part. 

(3) I^(e-|-*Fe): At finite temperature, this is replaced 
by I^;T(e + iFe) with 


T^.T{i+iTe) = P / dxF^{x,T)[l-F^{x,T)]Xx{eFiTe) ■ 

—00 

(B4) 

At |/i| — zq with 2zo = e -I- iFg for e > 0, 2^(e -I- fFe) 
diverges as (|/r| — zo)“^. Around e = 2|/r|, Re[I^(e -|- 
IFe)] has two minima ^ —F “^/2 around e ~ 2 |/r| ± -x/SFe 
and a maximum ~ T~^ at e = 2|/r|, while Im[I^(e -I- 
iFe)] has two extrema at e ~ 2 |/i| ± Te/'/i with values 
^ ±3-\/3/(8Fe). These indicate that this function varies 
very fast around e = 2|/r| for very small Fg. In a manner 
similar to the % function, at room temperature is 
a smooth function with respect to e > 0 for any Fg > 0 , 
and I^;T(*Fe) oc F^b 

In Fig. [ini (c) both I^(e -|- zFg) and X/i-rie + *Fg) are 
plotted for Fg = 10“^ eV and Fg = 0.03 eV at |^| = 
0.3 eV and T = 300 K. 
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